Visualizing hydrogren atom probability functions in 3D

Installing Plotly

import sys
!{sys.executable} -m pip install plotly

Requirement already satisfied: plotly in c:\users\molle\anaconda3\lib\site-packages (5.24.1)
Requirement already satisfied: tenacity>=6.2.0 in c:\users\molle\anaconda3\lib\site-packages (from plotly) (8.2.3)
Requirement already satisfied: packaging in c:\users\molle\anaconda3\lib\site-packages (from plotly) (25.0)

Importing libraries, defining constants

import numpy as np
import plotly.graph_objects as go
from plotly.subplots import make_subplots
from scipy import constants

a0 = constants.physical_constants['Bohr radius'][0]
Angstrom = constants.angstrom

Cartesian to spherical coordinates

def r_phi_theta(x, y, z):
    r = np.sqrt(x**2 + y**2 + z**2)
    phi = np.arctan2(y, x) + np.pi
    theta = np.arccos(z/r)
    return r, phi, theta

X, Y, Z = np.mgrid[-20*a0:20*a0:48j, -20*a0:20*a0:48j, -20*a0:20*a0:48j]
r_grid, phi_grid, theta_grid = r_phi_theta(X, Y, Z)

Normalized probability function

def probfunc(r, wavefunc):
    p = r**2 * np.conj(wavefunc) * wavefunc
    return p.real / np.max(p.real)

Question:

How come the wavefunctions n1l0m0 and n2l0m0 do not depend on $\theta$ and $\phi$?

Answer:

The wavefunctions n1l0m0 and n2l0m0 do not depend on $\theta$ and $\phi$ because the wavefunctions have no contribution from l or m (both values are zero). n is the principle quantum number, and is related only to the portion of the wavefunction dependent on the radius: R(r). l is the orbital angular momentum quantum number, which is related to both r and theta: R(r) and f($\theta$). m is the magnetic quantum number and is only related to the variable $\phi$: g($\phi$). Therefore, without the contributions of l and m in the wavefunction, it will not depend on $\theta$ and $\phi$, only r.

---

Question:

Which of the two probability functions plotted above is more centrally condensed, that is, a higher probability of the electron being close to the proton?

Answer:

Based on the two figures, the probability function for n1l0m0 points to being more centrally condensed, as it shows a sphere of smaller radius, meaning all of the electrons are in a fairly close orbit. This is opposed to the figure for n2l0m0, which shows the electrons being significantly more spread out around the nucleus, leading me to believe it is less centrally condensed.

---

Question:

Does the slice shown on the n1l0m0 plot match your expectations of the probability function?

Answer:

The slice shown on the n1l0m0 plot does match my expectations for the probability function. It demonstrates the spherical symmetry of the system, only changes as the radius changes (not with angles), and is a fairly simple function.

---

Question:

Based on what you see in the slices, how are the probability functions for n,l,m = 1,0,0 and n,l,m = 2,0,0 different?

Answer:

It looks like n1l0m0 continually decreases in probability the further out from the nucleus you get. On the other hand, n2l0m0's probability fluctates as it changes with the radius, and is not as straightforward of a system and probability function. Both form spheres, but n2l0m0 has a larger range and has more fluctuations in probability density.

---

Making a 3D isosurface rendering of n3l2m2 and n3l2m-2

Wavefunction

n3l2mneg2 = lambda r, phi, theta: (1/(162*np.sqrt(np.pi)*a0**1.5)) * (r**2/(a0**2)) * \
                               (np.exp(-r/(3*a0))) * (np.sin(theta)*np.sin(theta)) * \
                               (np.exp(-2j*phi))

Probability Function

prob_n3l2mneg2 = probfunc(r_grid, n3l2mneg2(r_grid, phi_grid, theta_grid))

n=3, l=2, m=-2, no slice

fig = go.Figure(data=go.Isosurface(
    x=X.flatten() / Angstrom, # x values of the grid in Angstroms
    y=Y.flatten() / Angstrom, # y values of the grid in Angstroms
    z=Z.flatten() / Angstrom, # z values of the grid in Angstroms
    value=prob_n3l2mneg2.flatten(), # independent variable
    isomin=0.05, # Minimum normalized probability density to render in an isosurface
    isomax=0.95, # Maximum normalized probability density to render in an isosurface
    opacity=0.4, # Set a low opacity so each surface is partially transparent
    colorscale='Plotly3_r', # Nice-looking color table
    surface_count=8, # number of isosurfaces to plot (2 by default: only min and max)
    colorbar_nticks=8, # colorbar ticks correspond to isosurface values
    caps=dict(x_show=False, y_show=False, z_show=False),
    #slices_x=dict(show=True, locations=[0]),
    ))

# Change axis lables and make the plot larger than the default
fig.update_layout(scene = dict(
                    xaxis_title='x (Angstroms)',
                    yaxis_title='y (Angstroms)',
                    zaxis_title='z (Angstroms)'),
                    width=700,
                    margin=dict(r=10, b=10, l=10, t=10))

fig.show()

n=3, l=2, m=-2, sliced

fig = go.Figure(data=go.Isosurface(
    x=X.flatten() / Angstrom, # x values of the grid in Angstroms
    y=Y.flatten() / Angstrom, # y values of the grid in Angstroms
    z=Z.flatten() / Angstrom, # z values of the grid in Angstroms
    value=prob_n3l2mneg2.flatten(), # independent variable
    isomin=0.05, # Minimum normalized probability density to render in an isosurface
    isomax=0.95, # Maximum normalized probability density to render in an isosurface
    opacity=0.4, # Set a low opacity so each surface is partially transparent
    colorscale='Plotly3_r', # Nice-looking color table
    surface_count=8, # number of isosurfaces to plot (2 by default: only min and max)
    colorbar_nticks=8, # colorbar ticks correspond to isosurface values
    caps=dict(x_show=False, y_show=False, z_show=False),
    slices_x=dict(show=True, locations=[0]),
    ))

# Change axis lables and make the plot larger than the default
fig.update_layout(scene = dict(
                    xaxis_title='x (Angstroms)',
                    yaxis_title='y (Angstroms)',
                    zaxis_title='z (Angstroms)'),
                    width=700,
                    margin=dict(r=10, b=10, l=10, t=10))

fig.show()